👉 Os testes de hipóteses em um contexto multivariado são muito mais complexos do que num contexto univariado.
Considerando a distribuição normal, temos
no caso univariado, uma média (\(\mu\)) e uma variância (\(\sigma^2\)), totalizando dois parâmetros,
no caso \(p\)-variado, \(p\) médias, \(p\) variâncias e \(\left( \begin{array}{c} p \\ 2 \end{array} \right)\) covariâncias, totalizando
\[p + p +\left( \begin{array}{c} p \\ 2 \end{array} \right) = \dfrac{1}{2} p(p+3)\]
parâmetros.
Por exemplo, para \(p = 10\), temos
\[\dfrac{1}{2} p(p+3) = \dfrac{10(10+3)}{2} = \dfrac{130}{2} = 65\]
parâmetros, para cada um dos quais uma hipótese poderia ser formulada.
Além disso, podemos estar interessados em testar hipóteses sobre subconjuntos desses parâmetros, ou ainda, sobre funções deles.
🤔 Por que utilizar testes multivariados ao invés de testes univariados?
\[\begin{eqnarray*} P(\text{pelo menos uma rejeição}) &=& 1 - P(\text{todos os 10 testes não rejeitarem } H_0) \\ &=& 1 - (0,95)^{10} = 0,40 \end{eqnarray*}\]
👉 Se as variáveis não forem independentes, teríamos \(0,05 < \alpha < 0,40\)
🤔 Por que utilizar testes multivariados ao invés de testes univariados?
Num primeiro momento, vamos revisar o caso univariado de se testar a média, \(\mu\), de uma população com \(\sigma^2\) conhecido.
Sejam as seguintes hipóteses:
\[H_0: \mu = \mu_0 \text{ vs. } H_a: \mu \neq \mu_0\]
Considere uma amostra aleatória \(X\) de tamanho \(n\), de forma que \(X \sim N(\mu, \sigma^2)\), com \(\sigma^2\) conhecido. Uma estatística de teste apropriada para este tipo de situação é
\[z = \dfrac{\bar{x}- \mu_0}{\dfrac{\sigma}{\sqrt{n}}} \stackrel{H_0}{\sim} N(0,1)\]
Em que \(\bar{x} = \dfrac{\displaystyle{\sum_{i=1}^n x_i}}{n}\). Para um nível \(\alpha\) de significância, rejeitamos \(H_0\) se \(|z| \geq z_{\frac{\alpha}{2}}\).
Equivalentemente, podemos usar \(z^2\) que segue uma distribuição \(\chi^2_1\) e rejeitamos a hipótese nula \(H_0\) se \(z^2 \geq \chi^2_1(\alpha)\).
Se \(n\) for grande, o Teorema Central do Limite garante que \(z\) é aproximadamente normal, mesmo que as observações não sejam normais.
Considere agora uma amoostra aleatória \(\mathbf{x} = \{X_1, X_2, \cdots, X_n\}\) de uma população \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\), com vetor de médias \(\mathbf{\mu} = \{\mu_1, \mu_2, \cdots, \mu_p\}^t\) e matriz de variâncias e covariâncias
\[\mathbf{\Sigma} = \begin{bmatrix} \sigma_{11} & \sigma_{12} & \cdots & \sigma_{1p}\\ \sigma_{21} & \sigma_{22} & \cdots & \sigma_{2p}\\ \vdots & \vdots & \ddots & \vdots\\ \sigma_{p1} & \sigma_{p2} & \cdots & \sigma_{pp} \end{bmatrix}\]
Queremos testar as seguintes hipóteses:
\[H_0: \mathbf{\mu} = \mathbf{\mu}_0 \text{ vs. } H_a: \mathbf{\mu} \neq \mathbf{\mu}_0\]
Podemos generalizar a estatistica \(z^2\) do teste univariado da seguinte forma:
\[Z^2 = n(\mathbf{\bar{x}}-\mathbf{\mu}_0)^t\mathbf{\Sigma}^{-1}(\mathbf{\bar{x}}-\mathbf{\mu}_0)\]
Se \(H_0\) é verdadeira, \(Z^2 \stackrel{H_0}{\sim} \chi^2_p\) e, portanto, rejeitamos a hipótese nula \(H_0\) se \(Z^2 > \chi^2_p(\alpha)\)
Exemplo 01: Considere uma amostra aleatória de tamanho \(n = 20\), com vetor de médias \(\mathbf{\bar{x}} = (71,45;164,7)^t\), de uma população normal bivariada cuja matriz de variâncias e covariâncias seja dada por
\[\mathbf{\Sigma} = \begin{bmatrix} 20 & 100 \\ 100 & 1000\\ \end{bmatrix}\]
Suponha que desejássemos testar as seguintes hipóteses
\[H_0: \mathbf{\mu} = [70, 170]^t \text{ vs. } H_a: \mathbf{\mu} \neq [70, 170]^t\]
Temos que,
\[\begin{eqnarray*} Z^2 &=& n(\mathbf{\bar{x}}-\mathbf{\mu}_0)^t\mathbf{\Sigma}^{-1}(\mathbf{\bar{x}}-\mathbf{\mu}_0) \\ &=& 20 \times \left[ \begin{array}{c} 71,45 - 70 & 164,7 - 170 \end{array} \right]\begin{bmatrix} 20 & 100 \\ 100 & 1000\\ \end{bmatrix}^{-1}\left[ \begin{array}{c} 71,45 - 70 \\ 164,7 - 170 \end{array} \right] \\ &=& 20 \times \left[ \begin{array}{c} 1,45 & -5,30 \end{array} \right]\begin{bmatrix} 0,1 & -0,01 \\ -0,01 & 0,002\\ \end{bmatrix}\left[ \begin{array}{c} 1,45 \\ -5,30 \end{array} \right] = 8,4026 \end{eqnarray*} \]
Usando um nível de significância \(\alpha = 0,05\), temos que \(\chi^2_{(2;0,05)} = 5,991\) e, portanto, rejeitamos a hipótese \(H_0: \mathbf{\mu} = [70, 170]^t\), uma vez que \(Z^2 = 8,4026 > 5,991\).
n = 20
p = 2
alpha = 0.05
xbarra = c(71.45,164.7)
mi0 = c(70,170)
Sigma = matrix(c(20,100,100,1000), nrow = p, byrow = TRUE)
# H0: mu = [70, 170]^t - Estatística Z2
Z2 = n*t(xbarra-mi0)%*%solve(Sigma)%*%(xbarra-mi0)
# valor crítico da distribuição qui-quadrado
qui = qchisq(1-alpha, p)
# p-valor
pvalor = 1-pchisq(Z2, p)
Z2 > qui # Rejeita H0Considere o caso univariado de uma amostra aleatória de tamanho \(n\), \(X_1, X_2, \cdots, X_n\), da distribuição \(N(\mu,\sigma^2)\), em que ambos os parâmetros \(\mu\) e \(\sigma\) são desconhecidos.
Considere também, que desejamos testar se um determinado valor \(\mu_0\) é um possível valor para a média populacional \(\mu\). Do ponto de vista de teste de hipóteses, isso é equivalente a testar:
\[H_0: \mu = \mu_0 \hspace{0.5cm} \text{ vs } \hspace{0.5cm} H_a: \mu \neq \mu_0\]
Neste caso, uma estatística de teste apropriada é
\[t = \displaystyle{\frac{\bar{X} - \mu_0}{s/\sqrt{n}}}\]
em que \(\bar{X} = \dfrac{\displaystyle{\sum_{j=1}^{n}} X_{j}}{n}\) e \(s^{2} = \dfrac{\displaystyle{\sum_{j=1}^{n}}(X_{j}-\bar{X})^2}{n-1}\).
Esta estatística de teste tem uma distribuição t-Student com \(n − 1\) graus de liberdade. Rejeita-se \(H_0\), se o valor observado de \(|t|\) exceder um valor específico da distribuição t-Student com \(n − 1\) graus de liberdade.
Rejeitar \(H_{0}\) quando \(|t|\) é grande é equivalente a rejeitar \(H_{0}\) quando o valor de \[\begin{equation} t^{2}=\frac{(\bar{X}-\mu_{0})^{2}}{s^2/{n}}= n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0}) \label{eq1} \end{equation}\] for grande.
Uma vez que \(\bar{X}\) e \(s^{2}\) são observados, rejeitamos \(H_{0}\) em favor a \(H_{a}\) a um nível \(\alpha\) de significância, se: \[n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0})>t^{2}_{n-1,{\alpha}/{2}}\]
Ou, de forma equivalente,
\[n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0})>F_{1, n-1}(\alpha)\]
Se \(H_{0}\) não for rejeitada, concluímos que \(\mu_{0}\) é um possível valor para a média populacional \(\mu\).
Considere agora o problema de se determinar se um dado vetor \(\mathbf{\mu}_{0}\), \(p \times 1\), é um possível valor para o vetor de médias \(\mathbf{\mu}\) de uma distribuição normal \(p\)-variada.
Uma generalização natural da distância quadrática encontrada para o caso multivariado é dada por:
\[\begin{equation}T^{2}=(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\left ( \frac{\mathbf{S}}{n} \right )^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})=n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})\end{equation}\]
em que
\[\bar{\mathbf{x}} = \frac{1}{n} \sum_{j=1}^{n}{\mathbf{X}_{j}}\]
\[\mathbf{S} = \frac{1}{n-1} \sum_{j=1}^{n}(\mathbf{X}_{j} -\bar{\mathbf{x}})(\mathbf{X}_j-\bar{\mathbf{x}})^{t}\]
\[\mathbf{\mu}_{0} = \left[ \begin{array}{c} \mu_{10}\\ \vdots \\ \mu_{p0}\\ \end{array} \right] \]
A estatística \(T^{2}\) é chamada de \(T^{2}\) de Hotelling.
Se o valor observado da estatística \(T^{2}\) é muito grande, isto é, se \(\bar{\mathbf{x}}\) está “muito longe” de \({\mathbf{\mu}}_{0}\), a hipótese \(H_{0}:{\mathbf{\mu}}={\mathbf{\mu}}_{0}\) é rejeitada.
Proposição 01: Sejam \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) uma amostra aleatória de uma distribuição \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\) e \(T^2 = n(\bar{\mathbf{x}}-\mathbf{\mu}_{0})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-\mathbf{\mu}_{0})\). Então,
\[\displaystyle{\frac{n - p}{(n - 1)p}}T^2 \stackrel{H_0}{\sim} F_{p,n-p}\]
Dessa forma, um teste para as hipóteses
\[H_{0}:\mathbf{\mu} = \mathbf{\mu}_{0} \text{ vs. } H_{a} :\mathbf{\mu} \neq \mathbf{\mu}_{0}\]
pode ser formulado.
Ao nível \(\alpha\) de significância, rejeita-se \(H_{0}\) em favor de \(H_{a}\) se
\[T^{2}=n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}}) > \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
de acordo com a proposição 01. \(F_{p,n-p} (\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p; (n−p)]\).
Exemplo 02: Suponha que a matriz de dados para uma amostra aleatória de tamanho \(n = 3\) de uma população normal bivariada seja dada por
\[\mathbf{X} = \begin{bmatrix} 9 & 9 \\ 10 & 6\\ 8 & 6 \end{bmatrix}\]
Avalie o valor observado de \(T^2\) para \(\mathbf{\mu}_0 = [9,5]^t\). Qual a distribuição amostral de \(T^2\) neste caso?
Temos que
\[\mathbf{\bar{x}} = \begin{bmatrix} \bar{x}_1 \\ \bar{x}_2 \end{bmatrix} = \begin{bmatrix} \dfrac{6 + 10 + 8}{3} \\ \dfrac{9 + 6 + 3}{3} \end{bmatrix} = \begin{bmatrix} 8 \\ 6 \end{bmatrix}\] Além disso,
\[s_{11} = \dfrac{(6-8)^2 + (10-8)^2 + (8-8)^2}{2} = 4\]
\[s_{12} = \dfrac{(6-8)(9-6) + (10-8)(6-6) + (8-8)(3-6)}{2} = -3\]
\[s_{22} = \dfrac{(9-6)^2 + (6-6)^2 + (3-6)^2}{2} = 9\]
De forma que,
\[\mathbf{S} = \begin{bmatrix} 4 & -3 \\ -3 & 9 \end{bmatrix}\]
e,
\[\mathbf{S}^{-1} = \dfrac{1}{(4 \times 9) - ((-3) \times (-3))}\begin{bmatrix} 9 & 3 \\ 3 & 4 \end{bmatrix} = \begin{bmatrix} \frac{1}{3} & \frac{1}{9} \\ \frac{1}{9} & \frac{4}{27} \end{bmatrix}\]
Sendo a estatística \(T^2\) de Hotelling dada por
\[\begin{eqnarray*}T^2 &=& n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}}) \\ &=& 3 \times \begin{bmatrix} 8-9 & 6-5 \end{bmatrix}\begin{bmatrix} \frac{1}{3} & \frac{1}{9} \\ \frac{1}{9} & \frac{4}{27} \end{bmatrix}\begin{bmatrix} 8-9 \\ 6-5 \end{bmatrix}\\ &=& 3 \times \begin{bmatrix} -1 & 1 \end{bmatrix}\begin{bmatrix} -\frac{2}{9} \\ \frac{1}{27} \end{bmatrix} = \dfrac{7}{9}\end{eqnarray*}\]
Antes da amostra ser selecionada, a distribuição de \(T^2\) é dada por
\[\dfrac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) = \frac{(3-1) \times 2}{(3-2)}F_{2,3-2} = 4F_{2,1}\]
Sendo assim, para testarmos as hipóteses
\[H_{0}:\mathbf{\mu} = \begin{bmatrix} 9 \\ 5 \end{bmatrix} \text{ vs. } H_{a} :\mathbf{\mu} \neq \begin{bmatrix} 9 \\ 5 \end{bmatrix}\]
Basta comparar o valor observado de \(T^2\) com o valor crítico da distribuição \(F_{2,1}\), ao nível de significância de \(\alpha = 0,05\), dado por \(4F_{2,1}(0,05) = 4 \times 199,5 = 798\)
\[T^2 = \dfrac{7}{9} = 0,778 << 798 = 4F_{2,1}(0,05)\] Logo, não existem evidências significativas a favor da rejeição de \(H_0\), ao nível de 5% de significância.
Exemplo 03: A transpiração de 20 mulheres saudáveis foi analisada. As observações são trivariadas, a saber: \(X_1\) - taxa de suor, \(X_2\) - conteúdo de sódio e \(X_3\) - conteúdo de potássio. Deseja-se testar, ao nível de significância de 10%, as seguintes hipóteses:
\[H_0: \mathbf{\mu} = \left[ \begin{array}{c} 4 \\ 50 \\ 10\\ \end{array} \right] \text{ vs } H_a: \mathbf{\mu} \neq \left[ \begin{array}{c} 4 \\ 50 \\ 10\\ \end{array} \right]\]
dados = read.table("https://raw.githubusercontent.com/tiagomartin/est022/refs/heads/main/dados/suor.dat", col.names = c("suor","sodio","potassio"))
dados %>% str()'data.frame': 20 obs. of 3 variables:
$ suor : num 3.7 5.7 3.8 3.2 3.1 4.6 2.4 7.2 6.7 5.4 ...
$ sodio : num 48.5 65.1 47.2 53.2 55.5 36.1 24.8 33.1 47.4 54.1 ...
$ potassio: num 9.3 8 10.9 12 9.7 7.9 14 7.6 8.5 11.3 ...
# Teste de normalidade multivariada
# H_0: dados seguem distribuição normal trivariada
mvshapiro_test(as.matrix(dados))
Generalized Shapiro-Wilk test for Multivariate Normality
data: as.matrix(dados)
MVW = 0.94446, p-value = 0.2567
p-valor = 0,2567 > 0,10 = \(\alpha\) → \(NRH_0\)
One Sample Hotelling T Square Test
Hotelling T Sqaure Statistic = 9.738773
F value = 2.905 , df1 = 3 , df2 = 17 , p-value: 0.0649
Descriptive Statistics
suor sodio potassio
N 20.00000 20.00000 20.000000
Means 4.64000 45.40000 9.965000
Sd 1.69687 14.13465 1.904641
Detection important variable(s)
Lower Upper Mu0 Important Variables?
suor 3.555292 5.724708 4 FALSE
sodio 36.364555 54.435445 50 FALSE
potassio 8.747476 11.182524 10 FALSE
\[\dfrac{n-p}{(n-1)p}T^2 = F_{p, n-p}\]
A rejeição da hipótese nula nos conduz ao problema de identificar quais variáveis são responsáveis pela rejeição dessa hipótese.
A região de confiança (RC) é uma generalização dos intervalos de confiança univariados.
A definição dessa região nos auxilia a identificar os componentes do vetor de parâmetros que são responsáveis pela rejeição de \(H_0\).
Seja \(\mathbf{\theta}\) um vetor de parâmetros populacional desconhecidos e \(\mathbf{\Theta}\) o conjunto de todos os possíveis valores de \(\mathbf{\theta}\). A RC é uma região de possíveis valores de \(\mathbf{\theta}\).
A RC é chamada de região com \(100(1 - \alpha)\%\) de confiança se, antes da amostra ser selecionada,
\[P(\mathbf{\theta} \in RC) = 1 - \alpha\]
Esta probabilidade é calculada sob o verdadeiro, porém desconhecido, valor de \(\mathbf{\theta}\).
A RC para o vetor de médias populacional \(\mathbf{\mu}\), quando se dispõe de uma amostra aleatória da distribuição \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\) é dada por,
\[P\left[ n(\bar{\mathbf{x}}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{X}}-{\mathbf{\mu}}) \leqslant \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) \right] = 1 - \alpha\]
de forma que,
\[RC = \left\lbrace \mathbf{\mu} \in \mathbf{R}^p | n(\bar{\mathbf{x}}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}}) \leqslant \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) \right\rbrace \]
em que \(F_{p,n-p} (\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p; (n−p)]\).
Temos que a região de confiança é dada pelo hiper elipsoide de eixos determinados pelos autovetores da matriz de covariância amostral \(\mathbf{S}\) e cujas medidas são proporcionais às raízes quadradas dos respectivos autovalores.
Para verificar se um dado vetor \(\mathbf{\mu}_0\) pertence à região de confiança, basta calcular o valor de
\[n(\bar{\mathbf{x}}-{\mathbf{\mu}_0})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_0})\]
e compará-lo com o valor de
\[\frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
Para \(p \geqslant 4\) não é possível representar visualmente a RC. No entanto, podemos calcular as medidas dos eixos do hiper elipsoide de confiança centrado em \(\mathbf{\bar x}\):
\[n(\mathbf{\bar x}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\mathbf{\bar x}-{\mathbf{\mu}}) \leqslant c^2 = \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
Os semi-eixos têm medida
\[\sqrt{\lambda_i} \displaystyle{\frac{c}{\sqrt{n}}} = \sqrt{\lambda_i} \sqrt{\displaystyle{\frac{(n-1)p}{n(n-p)}F_{p,n-p}(\alpha)}}\]
unidades na direção de \(\mathbf{e}_i\), em que \(\lambda_i\) e \(\mathbf{e}_i\) são os autovalores e autovetores de \(\mathbf{S}\).
Começando no centro \(\mathbf{\bar x}\), os eixos do elipsoide de confiança são
\[\pm \left( \sqrt{\lambda_i} \sqrt{\displaystyle{\frac{(n-1)p}{n(n-p)}F_{p,n-p}(\alpha)}}\right) \mathbf{e}_i\]
em que \(\lambda_i\) e \(\mathbf{e}_i\) são os autovalores e autovetores de \(\mathbf{S}\), \(i = 1, \cdots, p\).
Exemplo 04: Uma empresa quer verificar se a média do desempenho dos funcionários em dois critérios-chave está de acordo com o padrão estabelecido pela diretoria.
Os critérios avaliados são:
A diretoria espera que a média populacional dos funcionários seja:
\[\mathbf{\mu}_0 = \left[ \begin{array}{c} 80 \\ 7,5 \end{array} \right]\]
dados = read.csv("https://raw.githubusercontent.com/tiagomartin/est022/refs/heads/main/dados/dados_desempenho.csv", row.names = 1)
dados %>% str()'data.frame': 20 obs. of 2 variables:
$ Produtividade: num 90.5 87.3 69.3 84.2 83.7 ...
$ Satisfacao : num 7.9 6.88 6.15 7.09 6.95 6.62 5.38 7.38 8.07 5.75 ...
One Sample Hotelling T Square Test
Hotelling T Sqaure Statistic = 67.55705
F value = 32.001 , df1 = 2 , df2 = 18 , p-value: 1.18e-06
Descriptive Statistics
Produtividade Satisfacao
N 20.000000 20.000000
Means 83.580500 6.431500
Sd 9.732401 1.026195
Detection important variable(s)
Lower Upper Mu0 Important Variables?
Produtividade 77.619031 89.541969 80.0 FALSE
Satisfacao 5.802916 7.060084 7.5 *TRUE*
p-valor = 0,00000118 < 0,05 = \(\alpha\) → \(RH_0\)
🤔 Qual variável foi a responsável?
Se a média esperada \(\mathbf{\mu}_0\) estiver fora da elipse, há indícios de que o desempenho real dos funcionários não está de acordo com a expectativa.
x_barra = colMeans(dados)
n = nrow(dados)
p = ncol(dados)
S = var(dados)/n
alpha = 0.05
c2 = (n-1)*p*qf(1-alpha,p,n-p)/(n-p)
ci.prod=ci_mean(dados$Produtividade, type = "t")$interval
ci.satis=ci_mean(dados$Satisfacao, type = "t")$interval
# Cálculo da elipse de confiança (nível 95%)
Radius = sqrt(eigen(S)$values) * sqrt(c2)
segments = 51
angles = (0:segments) * 2 * pi/segments
unit.circle = cbind(cos(angles), sin(angles))
Q = sweep(eigen(S)$vectors, 2, Radius, "*")
elipse_conf = data.frame(sweep(unit.circle %*% t(Q), 2, as.vector(x_barra), "+"))
names(elipse_conf) = c("Produtividade", "Satisfacao")
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept=ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Quando \(p \geqslant 4\), não conseguimos visualizar graficamente a região de confiança, de forma que necessitaremos de um procedimento alternativo para descobrir qual componente do vetor de parâmetros foi o responsável pela rejeição de \(H_0\).
A ideia básica é fazermos intervalos de confiança para cada componente, ou para combinações lineares dos componentes do vetor de parâmetros.
Sejam \(\mathbf{x} \sim N_p(\mathbf{\mu}, \mathbf{\Sigma})\) e a combinação linear
\[Z = \mathbf{l}^t \mathbf{x} = l_1X_1 + l_2X_2 + \cdots + l_pX_p\]
Temos que
\[\mu_Z = E(Z) = \mathbf{l}^t \mathbf{\mu} \hspace{0.5cm} \text{e} \hspace{0.5cm} \sigma^2_Z = \mathbf{l}^t \mathbf{\Sigma} \mathbf{l}\]
de forma que
\[Z \sim N(\mathbf{l}^t \mathbf{\mu}, \mathbf{l}^t \mathbf{\Sigma} \mathbf{l})\]
Logo, se \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) é uma amostra aleatória de uma população \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\), segue que \(Z_1, Z_2, \cdots, Z_n\), definidos por
\[Z_i = \mathbf{l}^t \mathbf{x}_i, \hspace{0.5cm} i = 1, 2, \cdots n\]
é uma amostra aleatória de uma população \(N(\mathbf{l}^t \mathbf{\mu}, \mathbf{l}^t \mathbf{\Sigma} \mathbf{l})\).
Da teoria normal univariada, para \(\mathbf{l}\) fixo e \(\sigma^2\) desconhecido, temos que um intervalo de \(100(1-\alpha)\%\) de confiança para \(\mu_Z = \mathbf{l}^t \mathbf{\mu}\) é baseado na razão t-Student, isto é,
\[ t = \displaystyle{\frac{\bar{z} - \mu_Z}{s_Z/\sqrt{n}}} = \displaystyle{\frac{\sqrt{n}(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})}{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}} \]
O que conduz a,
\[\bar{z} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{s_Z}}{\sqrt{n}}} \leqslant \mu_Z \leqslant \bar{z} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{s_Z}}{\sqrt{n}}}\]
ou,
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Podemos construir vários intervalos de confiança sobre combinações lineares dos componentes do vetor \(\mathbf{\mu}\), cada um associado com um coeficiente de confiança \(1 − \alpha\), escolhendo diferentes vetores de constantes \(\mathbf{l}\). Porém, o coeficiente de confiança conjunto do grupo de intervalos resultantes não seria mais \(1 − \alpha\).
É desejável então, associar um coeficiente de confiança coletivo de \(1 − \alpha\) aos intervalos de confiança que podem ser gerados para diferentes escolhas de \(\mathbf{l}\) (Bônus).
Ônus: intervalos que têm amplitudes maiores (mais largos, menos precisos) do que os intervalos apresentados anteriormente via a distribuição amostral t-Student com \(n − 1\) graus de liberdade.
Uma solução para o problema de cobertura é a construção dos intervalos simultâneos. Dado um conjunto de dados observados \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) e uma particular escolha do vetor \(\mathbf{l}\), temos que o intervalo de confiança
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
representa o conjunto de valores de \(\mathbf{l}^t \mathbf{\mu}\) para os quais
\[|t| = \left|\displaystyle{\frac{\sqrt{n}(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})}{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}}\right| \leqslant t_{\left[ (n-1);\frac{\alpha}{2}\right]}\]
Ou, equivalentemente,
\[t^2 = \displaystyle{\frac{n(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})^2}{\mathbf{l}^t \mathbf{S} \mathbf{l}}} = \displaystyle{\frac{n\left[ \mathbf{l}^t( \bar{\mathbf{x}} - \mathbf{\mu})\right] ^2}{\mathbf{l}^t \mathbf{S} \mathbf{l}}} \leqslant t^2_{\left[ (n-1);\frac{\alpha}{2}\right]}\]
Uma região de confiança simultânea é dada pelos valores de \(\mathbf{l}^t \mathbf{\mu}\) tais que \(t^2\) é relativamente pequeno para todas as escolhas de \(\mathbf{l}\).
Parece razoável esperar que a constante \(t_{\left[ (n-1);\frac{\alpha}{2}\right]}\) seja trocada por um valor maior, \(c^2\), quando as declarações forem desenvolvidas para muitas escolhas de \(\mathbf{l}\).
A constante \(c^{2}\) deve ser obtida para manter o valor global do coeficiente de confiança \(1-\alpha\) simultâneo para todas as escolhas de \(\mathbf{l}\).
Pode ser demonstrado que, tomando
\[c^{2}= \displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}\]
tem-se a garantia do valor global para o coeficiente de confiança igual a \(100(1-\alpha)\%\) para qualquer escolha não-nula do vetor \(\mathbf{l}\).
É conveniente nos referirmos aos intervalos simultâneos como intervalos \(T^2\), uma vez que a probabilidade de cobertura é determinada pela distribuição de \(T^2\)
Assim, os intervalos de confiança simultâneos (ou intervalos \(T^2\)) são dados por:
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - \sqrt{\displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}} \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + \sqrt{\displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}} \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Exemplo 05: Voltando ao exemplo anterior referente à uma empresa que quer verificar se a média do desempenho dos funcionários em dois critérios-chave (produtividade e satisfação) está de acordo com o padrão estabelecido pela diretoria.
Encontre os intervalos de confiança simultâneos \(T^2\) de Hotelling e compare com os intervalos univariados.
lim = cbind(x_barra - sqrt(c2)*sqrt(diag(S)),
x_barra + sqrt(c2)*sqrt(diag(S)))
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept = ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
geom_hline(yintercept = lim[2,], linetype="dashed", color = "blue") +
geom_vline(xintercept = lim[1,], linetype="dashed", color = "blue")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Uma alternativa aos intervalos simultâneos são os Intervalos de Confiança com Correção de Bonferroni, baseada na Desigualdade de Bonferroni.
Sejam \(A_1, A_2, \cdots, A_m\) uma coleção de eventos em um espaço de probabilidade tais que \(P(A_i) = 1 - \alpha_i\), \(i = 1, \cdots, m\). Então,
\[ \begin{eqnarray*} P(\cap_{i = 1}^m A_i) &=& 1 - P \left(``\text{pelo menos um } A_i \text{ não ocorre}"\right) \\&\geqslant& 1 - \displaystyle{\sum_{i = 1}^m P(A_i^c)} = 1 - \displaystyle{\sum_{i = 1}^m \alpha_i} \\&=& 1 - (\alpha_1 + \alpha_2 + \cdots + \alpha_m) \end{eqnarray*} \]
A desigualdade de Bonferroni nos permite o controle da taxa de erro global \(\alpha_1 + \alpha_2 + \cdots + \alpha_m\), independente da estrutura de correlação que está por trás dos intervalos de confiança.
Assim, se o problema envolve a construção de \(m\) intervalos de confiança, a ideia é fazer
\[\alpha^{\ast}=\dfrac{\alpha}{m}\]
e tomar os intervalos separadamente, dados por
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha^{\ast}}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha^{\ast}}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Observe agora, que o coeficiente de confiança coletivo vale pelo menos,
\[1 - \left( \underbrace{\displaystyle{\frac{\alpha}{m}} + \displaystyle{\frac{\alpha}{m}} + \cdots + \displaystyle{\frac{\alpha}{m}}}_{\text{m termos}}\right) = 1 - \alpha\]
Estes intervalos podem, de forma mais apropriada, serem comparados com os intervalos \(T^2\).
Exemplo 06: Voltando ao exemplo anterior referente à uma empresa que quer verificar se a média do desempenho dos funcionários em dois critérios-chave (produtividade e satisfação) está de acordo com o padrão estabelecido pela diretoria.
Encontre os intervalos de confiança com correção de Bonferroni e os compare com os intervalos simultâneos \(T^2\) de Hotelling e com os intervalos univariados.
m=2
alpha_corr = alpha/m
ci_bonf.prod=ci_mean(dados$Produtividade, type = "t", probs=c(alpha_corr/2, 1-alpha_corr/2))$interval
ci_bonf.satis=ci_mean(dados$Satisfacao, type = "t", probs=c(alpha_corr/2, 1-alpha_corr/2))$interval
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept = ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
geom_hline(yintercept = lim[2,], linetype="dashed", color = "blue") +
geom_vline(xintercept = lim[1,], linetype="dashed", color = "blue") +
geom_hline(yintercept = ci_bonf.satis, linetype="dashed", color = "#006400") +
geom_vline(xintercept = ci_bonf.prod,linetype="dashed", color = "#006400")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Novamente, vamos revisar o caso univariado. Sejam uma amostra aleatória \(x_{11}, x_{12}, \cdots, x_{1n_1}\) de uma população \(N(\mu_1, \sigma^2_1)\) e uma segunda amostra \(x_{21}, x_{22}, \cdots, x_{2n_2}\) de uma população \(N(\mu_2, \sigma^2_2)\). Vamos assumir que as amostras são independentes e que \(\sigma^2_1 = \sigma_2^2 = \sigma^2\), com \(\sigma^2\) desconhecido.
Sejam as seguintes hipóteses:
\[H_0: \mu_1 = \mu_2 \text{ vs. } H_a: \mu_1 \neq \mu_2\]
Uma estatística de teste apropriada para este tipo de situação é
\[t = \dfrac{\bar{x_1}- \bar{x}_2}{s_p\sqrt{\dfrac{1}{{n_1}}+\dfrac{1}{n_2}}} \stackrel{H_0}{\sim} t_{n_1 + n_2 - 2}\]
em que
\[s_p^2 = \dfrac{(n_1-1)s_1^2 + (n_2-1)s_2^2}{n_1 + n_2 - 2}\]
Considere agora uma amostra aleatória \(\mathbf{x_1} = \{X_{11}, X_{12}, \cdots, X_{1n_1}\}\) de uma população \(N_p(\mathbf{\mu_1}, \mathbf{\Sigma_1})\) e uma segunda amostra \(\mathbf{x_2} = \{X_{21}, X_{22}, \cdots, X_{2n_2}\}\) de uma população \(N_p(\mathbf{\mu_2}, \mathbf{\Sigma_2})\)
Assumiremos que as duas amostras são independentes e que \(\mathbf{\Sigma_1} = \mathbf{\Sigma_2} = \mathbf{\Sigma}\).
Queremos testar as seguintes hipóteses:
\[H_0: \mathbf{\mu}_1 = \mathbf{\mu}_2 \text{ vs. } H_a: \mathbf{\mu}_1 \neq \mathbf{\mu}_2\]
Quandoo as matrizes de covariâncias populacionais são iguais, as hipóteses podem ser testadas utilizando a seguinte estatística,
\[T^2 = (\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)^t \left[ \left( \dfrac{1}{n_1} + \dfrac{1}{n_2} \right) \mathbf{S}_p \right]^{-1}(\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)\]
em que
\[\mathbf{S}_p = \dfrac{(n_1-1)\mathbf{S}_1 + (n_2-1)\mathbf{S}_2}{n_1 + n_2 - 2}\]
Sob \(H_0\), temos que
\[T^2 \stackrel{H_0}{\sim} \displaystyle{\frac{p(n_1 + n_2 - 2)}{n_1 + n_2 - p - 1}} F_{(p,n_1+n_2-p-1)}\]
Ao nível \(\alpha\) de significância, rejeita-se \(H_{0}\) em favor de \(H_{a}\) se
\[T^{2}\geq \displaystyle{\frac{p(n_1 + n_2 - 2)}{n_1 + n_2 - p - 1}} F_{(p,n_1+n_2-p-1)}(\alpha)\]
em que \(F_{(p,n_1+n_2-p-1)}(\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p,(n_1+n_2-p-1)]\).
\[\dfrac{n_1+n_2-p-1}{(n_1+n_2-2)p}T^2 = F_{p, n_1+n_2-p-1}\]
Para o caso em que \(\mathbf{\Sigma}_1 \neq \mathbf{\Sigma}_2\), ou seja, as matrizes de covariâncias populacionais se diferem, sugere-se o uso da estatística
\[T^2 = (\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)^t \left[ \left( \dfrac{\mathbf{S}_1}{n_1} + \dfrac{\mathbf{S}_2}{n_2} \right) \right]^{-1}(\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)\]
Para grandes amostras, pode-se utilizar
\[T^2 \stackrel{H_0}{\sim} \chi^2_p,\] sendo \(\chi^2_p\) a distribuição qui-quadrado com \(p\) graus de liberdade.
A hipótese nula deve ser rejeitada se \(T^2 \geq \chi ^2_p(\alpha)\) para o nível de significância \(\alpha\).